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A deformed relativistic Hartree Bogoliubov (RHB) theory in continuum is developed aiming at 
a proper description of exotic nuclei, particularly those with a large spatial extension. In order 
to give an adequate consideration of both the contribution of the continuum and the large spatial 
distribution in exotic nuclei, the deformed RHB equations are solved in a Woods-Saxon (WS) basis in 
which the radial wave functions have a proper asymptotic behavior at large distance from the nuclear 
center. This is crucial for the proper description of a possible halo. The formalism of deformed 
RHB theory in continuum is presented. A stable nucleus, ??Mg and a weakly-bound nucleus, “*Mg, 
are taken as examples to present numerical details and to carry out necessary numerical checks. 
In addition, the ground state properties of even-even magnesium isotopes are investigated. The 
generic conditions of the formation of a halo in weakly bound deformed systems and the shape of 
the halo in deformed nuclei are discussed. We show that the existence and the deformation of a 
possible neutron halo depend essentially on the quantum numbers of the main components of the 
single particle orbitals in the vicinity of the Fermi surface. 


PACS numbers: 21.60.-n, 21.10.-k, 21.10.Gv, 21.60.Jz 


I. INTRODUCTION 


The development of radioactive ion beam facilities 
around the world [1-7] stimulates very much the study 
of nuclei far from the 8 stability line [8-16]. Some 
new and exotic phenomena have been observed in nuclei 
close to drip lines such as neutron or proton halos [17— 
19], changes of nuclear magic numbers [20], pygmy reso- 
nances [21], etc. In halo nuclei, the fact of an extremely 
weakly binding leads to many new features, e.g., the 
coupling between bound states and the continuum due 
to pairing correlations and very extended spatial density 
distributions. Therefore one must consider properly the 
asymptotic behavior of nuclear densities at large distance 
r from the center and treat in a self consistent way the 
discrete bound states, the continuum and the coupling 
between them in order to give a proper theoretical de- 
scription of the halo phenomenon [22-24]. This could 
be achieved by solving the non-relativistic Hartree-Fock- 
Bogoliubov (HFB) [25, 26] or the relativistic Hartree Bo- 
goliubov (RHB) [27-29] equations in coordinate (r) space 
which can fully take into account the mean-field effects 
of the coupling to the continuum. The resonant-BCS 
(rBCS) approach presents an other method to include 
the contribution of the resonant continuum which has 
been used to study halo phenomena [30, 31]. 

The solution of the coupled differential equations of 
HFB and RHB theories is particular simple in spheri- 
cal systems with local potentials, where one-dimensional 
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Numerov or Runge-Kutta methods [32] can be applied 
and this is true even for non-local problems where Finite 
Element Methods (FEM) [33] have been used. A differ- 
ent method to solve such equations is the expansion of 
the single particle wave functions in an appropriate ba- 
sis. The oscillator basis has been used for this purpose 
with a great success for deformed or non-local systems in 
the past [34-37]. The Woods-Saxon basis has been pro- 
posed in Ref. [38] as a reconciler between the harmonic 
oscillator basis and the integration in coordinate space. 
Woods-Saxon wave functions have a much more realistic 
asymptotic behavior at large r than the harmonic oscil- 
lator wave functions do. A discrete set of Woods-Saxon 
wave functions is obtained by using box boundary condi- 
tions to discretize the continuum. It has been shown in 
Ref. [38] for spherical systems that the solution of the rel- 
ativistic Hartree equations in a Woods-Saxon basis is al- 
most equivalent to the solution in coordinate space. The 
Woods-Saxon basis has also been used in more compli- 
cated situations, e.g., for the description of exotic nuclei 
where both deformation and pairing have to be taken 
into account. Recently, for spherical systems, both non- 
relativistic and relativistic Hartree-Fock-Bogoliubov the- 
ories with forces of finite range have been investigated in 
a Woods-Saxon basis [39, 40]. 

Over the past years, lots of efforts have been made to 
develop a deformed relativistic Hartree (RH) theory [41] 
and a deformed relativistic Hartree Bogoliubov theory in 
continuum [42]. As a first application, halo phenomena 
in deformed nuclei have been investigated within the con- 
tinuum RHB theory and some brief results can be found 
in Ref. [43]. In this paper we present the full version of 
the theoretical framework with all the details. 
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Spherical symmetry facilitates considerably the treat- 
ment of the continuum in non-relativistic HFB [25, 26] 
and in relativistic RHB theory [27-29] in r-space. Since 
most of the known nuclei are deformed, interesting ques- 
tions arise, whether or not deformed halos exist and what 
new features can be expected in deformed exotic nu- 
clei [44-49]. Such questions can be answered by the de- 
formed counterparts of the HFB or RHB theories in coor- 
dinate space. From the experimental point of view, ?4Ne 
is measured to be a strongly deformed halo nucleus [50], 
and for the well deformed magnesium isotopes, °°Mg is 
probably a halo nucleus too [51]. Nevertheless for de- 
formed nuclei, to solve the HFB or RHB equations in 
r space becomes much more sophisticated and numeri- 
cally very time consuming. Many efforts have been made 
to develop non-relativistic HFB theories either in (dis- 
cretized) coordinate space or in a scaled oscillator basis 
with improved asymptotic behavior [37]. The HFB equa- 
tions have been solved in three-dimensional coordinate 
space by combining the imaginary time approach and the 
two basis method [52] with a truncated basis composed of 
discrete localized states and discretized continuum states 
up to a few MeV [53]. Alternatively, the HFB equa- 
tions have been solved on a two-dimensional basis-spline 
Galerkin lattices [54-56] or on a three-dimensional Carte- 
sian mesh [57] using the canonical-basis approach [58]. 
Recently, the Gaussian expansion method is used to solve 
the HF and HFB equations for deformed nuclei [59] and 
continuum Skyrme-Hartree-Fock-Bogoliubov approaches 
have been developed both for spherical and deformed nu- 
clei [60]. The deformed relativistic Hartree Bogoliubov 
(RHB) theory has only been solved in the conventional 
harmonic oscillator basis [61-64] and neither the above- 
mentioned approaches nor other methods which could 
improve the asymptotic behavior of the nuclear densities 
at large r have been implemented in the deformed RHB 
theory so far. 


In this paper we present a method, which allows to 
take into account at the same time the coupling to the 
continuum, deformations, and pairing correlations in a 
fully self-consistent way. For this purpose we expand the 
deformed Dirac spinors in a basis of spherical Dirac wave 
functions obtained by the solution of the Dirac equations 
for potentials with spherical Woods-Saxon shape. This 
idea is similar to a method proposed in Ref. [65] for the 
solution of the deformed relativistic mean field (RMF) 
equations in light nuclei, where the deformed Dirac- 
spinors were expanded in terms of the self-consistent so- 
lutions of the spherical RMF-equations. As compared to 
these early calculations our method is simpler, because it 
is based on Woods-Saxon wave functions. On the other 
side it is more general, because it allows to include pairing 
correlations, which play an essential role in the formation 
of halo structures. 


The paper is organized as follows. In Sec. II, we give 
the formalism of the deformed RHB theory in continuum. 
The numerical details are presented in Sec. III and we 
discuss applications and detailed results for magnesium 


isotopes in Sec. IV. A summary is given in Sec. V. 


II. FORMALISM OF THE DEFORMED 
RELATIVISTIC HARTREE BOGOLIUBOV 
THEORY IN CONTINUUM 


The starting point of relativistic mean field theory is 
a Lagrangian density where nucleons are described as 
Dirac spinors which interact via the exchanges of effective 
mesons (a, w, and p) and the photon [66-72], 


L=4(iðĝ— M) Y+ 5.08! — U (0) — gobo 


1 1 - 
~ ZW + mowa — gubo 


2 
lz PHV Las zu 7 F 
~ gE R” + 5m Pu — goh prh 
1 D =1= T3 
= qe ~~ ep 9 Ar, (1) 


where M is the nucleon mass, and Mo, Jo, Mw, Jus Mp, Ip 
masses and coupling constants of the respective mesons. 
The nonlinear self-coupling for the scalar meson is given 
by [73] 

U(o) = imzo? + 2a + Pot, (2) 
and field tensors for the vector mesons and the photon 
fields are defined as 

Quy = OMEN = rwu, 
Ruy = Ou Pv = vpu = Ip (Pu x pv); (3) 

Fy = Or Ay = Anz. 

Pairing correlations are crucial in the description of 
open shell nuclei. For exotic nuclei, the conventional 
BCS approach turns out to be only a poor approxima- 
tion [26]. Starting from the Lagrangian density (1), a rel- 
ativistic theory of pairing correlations in nuclei has been 
developed by Kucharek and Ring [74]. If we neglect the 
Fock terms as it is usually done in the covariant density 
functional theory, the Dirac Hartree Bogoliubov (RHB) 
equation for the nucleons reads, 


few (ea py) (0) ame (%) a0 


where Ex is the quasiparticle energy, A is the chemical 
potential, and hp is the Dirac Hamiltonian, 


hp(r,r')=a-p+V(r)+hM+S(r)). (5) 
The scalar and vector potentials 


S(r) = goo(r), (6) 


LEA), (7) 


depend on the scalar field o and on the time-like compo- 
nents w, p°, and A? of the iso-scalar vector field w, the 


V(r) = gow (r) + gpT3p°(r) +e 
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3-component of iso-vector vector field p and the photon 
field. 
The equations of motion for the mesons and the photon 


(-A +8,U(a)) o(r) = —gops(n), 
-A +m?) w (r) = gupr(r), (8) 
(~A +m) lr) = gpp3(r), 
~KA%"r) = epr), 


have as sources the various densities 
ps(r) = X Vi (r)yoVe(r), 
k>0 
polr) = E Vi Vr), 
k>0 


plr) = D Vilr)rsVelr), ©) 
pelr) = D Vivl), 
k>0 


where, according to the no-sea approximation, the sum 

over k > 0 runs over the quasi-particle states correspond- 

ing to single particle energies in and above the Fermi sea. 
The pairing potential reads, 


s3p2 

: (en a ee 

A(r1s1p1,7282p2) = X VPP (T1, T2; $1P1, S2P2, S1P1, S2P2) 
f a 
$ıPı 


X K(T1S1P1, T282P2), (10) 


where p = 1,2 is used to represent the large and small 
components of the Dirac spinors. VPP is the effective 
pairing interaction and «(r1s1p4, r2s5%p3) is the pairing 
tensor [75]. 

In the particle-particle (pp) channel, we use a density 
dependent zero range force, 


Pipens Vox (1—P?)6(e1—r2) (1 _ A (1) 


4(1 — P?) projects onto spin S = 0 component in the 
pairing field. In this case the gap equation (10) has the 
simple form 


A(r) = Vo(1 = p(r)/psat)e(r), 


and we need only the local part of the pairing tensor 


n(r) = $ Vi (r)U,(r), 


k>0 


(12) 


(13) 


Details of the calculation of the pairing interaction and 
the pairing tensor are given in Appendices B and E re- 
spectively. 

For axially deformed nuclei with spatial reflection 
symmetry, we expand the potentials S(r) and V(r) in 
Eqs. (6) and (7) and the densities in Eq. (9) in terms of 
the Legendre polynomials [76], 


f(r) = XL fal) Pa (cos0), A =0,2,4, +, (14) 
A 


with 


2A\4+1 
f(r) === 


J df (r) P(Q). (15) 


The quasiparticle wave functions U;, and Vp in Eq. (4) 
are Dirac spinors. Each of them is expanded in terms 
of spherical Dirac spinors Ynsm(rsp) with the eigenval- 
ues Eng Obtained from the solution of a Dirac equation 
ne) containing spherical potentials S% (r) and V(r) 
of Woods-Saxon shape [38, 77]: 


U;.(rsp) = 5 U Oey Pram (PSP); (16) 
Vi(rsp) = Y ve?) Prem (TSP). (17) 
The basis wave function reads 
1 iGna(T)Yjim (Qs) 
ngem = 7 ) 1 
p (rs) m ( —Farn(r)Y ln (Qs) ( 8) 


where G,,(r)/r and F,,,(r)/r the radial wave func- 
tions for the upper and lower components. The spher- 
ical spinor Ynem is characterized by the radial quan- 
tum number n, angular quantum j and the parity m. j 
and m are combined to the relativistic quantum number 
k = m(—1)/+1/2(j + 1/2) which runs over positive and 
é Yin and Yon are the 
spinor spherical harmonics where | = j + $sign(«) and 
l= j- ssign(s). 

Pnem(rsp) is the time reversal state of Yrem(Trsp). 
These states form a complete spherical and discrete basis 
in Dirac space (see Appendix A for details). Because of 
the axial symmetry the z-component m of the angular 
momentum j is a conserved quantum number and the 
RHB Hamiltonian can be decomposed into blocks char- 
acterized by m and parity a. For each mz-block, solving 
the RHB equation (4) is equivalent to the diagonalization 
of the matrix 


negative integers A = +1, +2,- 


A-X B Un\ _ Uk 
TENORON. 
where 

Ue = (uho) Ve = (kew) 20) 
and 


A= Cee a ((nam|hp|n'r’, m)) ’ (21) 


B= Ce 


Cee) = ((nam|A|n’k’,m)) . (22) 
Further details are given in Appendix B. 

Since we use a zero range pairing force we have to 
introduce a pairing cutoff in the sums of Eqs. (9) and 


(13) over the quasiparticle space. In the present work, a 
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smooth cut off is adopted where two parameters, Eo)" 


and 24, are introduced and the square root of the factor 


1 Ey, — Et? 
s(Ex) = 3 | 1- Sas] > 3) 
2 (Ex =H ) T (Tiit ) 


cut 


is multiplied in the occupation component V(r) of each 

quasi particle state with v? < 1/2. Note that this smooth 

cutoff is similar as the soft cutoff proposed in Ref. [78]. 
The total energy of a nucleus is 


E = Enueleon } Eg H Ew } Ep } Be H Bem 
= XL 2A — Bx) ug — Epair 
k 


-i fèr [goo (r)ps(r) + U(o)] 
-5 J Prguw(r)ps(r) 


-5 / d’rg,p(r) [ec (r) = py (r)] 


1 
= / Pr Agp? (r) + Eom. (24) 
where 
2 
v = Jervie) = 5 a) . (25) 


For a zero range force the pairing field A(r) is local 
and the pairing energy is calculated as 


1 
Epas = =F | Pre(r) AC). (26) 
The center of mass correction energy 
1 
Ecm. = — (P? 2 
m= yo (P*), (27) 


is calculated after variation with the wave functions of 
the self-consistent solution [79, 80] or in the oscillator 
approximation 


Eom. = = x 41 x Ae MeV, (28) 


Details are given in Appendix G. The root mean square 
(rms) radius is calculated as 


= (/ dr retro!) (29) 


where T represents the proton, the neutron, or the nu- 
cleon. The rms charge radius is calculated simply as 
ra, = r2 + 0.64 fm?. The intrinsic multipole moment 
is calculated by 


Qr,A = 
a fo [rto a(r]. (30) 


The quadrupole deformation parameter is obtained from 
the quadrupole moment by 


(31) 


where N, refers to the number of neutron, proton, or 
nucleon. 


Ill. NUMERICAL DETAILS AND ROUTINE 
CHECKS 


A. Details on the Woods-Saxon basis 


For numerical reasons several parameters have to be 
introduced in the calculations, e.g., the mesh size Ar, 
the box size Rpox for the determination of the basis wave 
functions by solving the spherical Dirac equations with 
the Hamiltonian hO, the maximal A-value Amax in the 
expansion Eq. (14) of the deformed fields and densities, 
the cutoff parameters for the radial and angular quantum 
numbers n and « in the expansion of Eqs. (16) and (17), 
Nmax and Kmax. Instead of Nmax, we introduced an energy 
cutoff parameter Et,, for positive energy states in the 
Woods-Saxon basis and in each «-block, the number of 
negative energy states in the Dirac sea is the same as 
that of positive energy states above the Dirac gap [38]. 
We have investigated the dependence of our results on 
these parameters in spherical and deformed relativistic 
Hartree models [38, 41]. It is found that a box of the size 
Rpox = 4ro A13 with ro = 1.2 fm, a step size Ar = 0.1 
fm, Amax = 4, and |Kmax| = 15 leads in light nuclei to 
an acceptable accuracy of less than 0.1 % for the binding 
energies, the rms radii, and the quadrupole moments . 

In the present work we use the determination of the 
Woods-Saxon basis a box size Rpox = 20 fm, a mesh size 
Ar = 0.1 fm and a cutoff energy Et, = 100 MeV. In 
each «-block in the Woods-Saxon basis, the number of 
negative energy states in the Dirac sea is the same as 
that of positive energy states above the Dirac gap. In 
Sec. IIIC we investigate the convergence of our results 
with respect to these three parameters. 

In order to reduce the computational time, Amax = 4 
and |Kmax| = 10 are used in this work. The parameter 
sets NL3 [81] and PK1 [79] are used for the Lagrangian 
density. Note that the center of mass correction energy is 
calculated differently with these two parameter sets. For 
NL3, the empirical formula in Eq. (28) is used and for 
PK1, the center of mass correction energy is calculated 
microscopically (see Appendix G). 


B. Parameters for the pairing force 


There are two parameters in the phenomenological 
pairing force Eq. (11), namely, Vo and psat, and two ad- 
ditional ones in the smooth cutoff Eq. (23). We take 
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TABLE I. Determination of the parameters for the pairing 
force used in the deformed RHB calculations presented in this 
work. In the last column is given the proton pairing energy 
E5air from the SRHBHO and deformed RHB calculations for 
the spherical nucleus 7°Mg. 


Model 


SRHBHO DIS [83] 


Surface 6 


RHB 
with Psat = 0.152 fm~? 
smooth cutoff] PLR = 60 MeV 


Tt? = 5.65 MeV 


cut 


TABLE II. Ground state properties of 7??Mg from deformed 
RHB calculations with different cutoff parameters in the 
Woods-Saxon basis compared with the results of spherical 
RCHB [29] theory. 


deformed RHB 
100 200 
—0.8992 —0.9072 
2.3823 2.3866 
2.5902 
3.0052 
—9.2294 


—0.9063 
2.3871 
2.5900 
3.0049 

—9.2381 


2.5910 
3.0073 
—9.1165 


the empirical value 0.152 fm™? for the saturation den- 
sity Psat. The pairing strength Vo, the cutoff EX? is 
fixed by reproducing the proton pairing energy of the 
Gogny force D1S in the spherical nucleus 7??Mg. We 
first calculate the ground state properties of ?9Mg by 
using the spherical relativistic Hartree Bogoliubov the- 
ory in a harmonic oscillator basis (SRHBHO) [82] in 
which the Gogny-D1S [83] force is used in the pp channel. 
The pairing energy for protons is obtained as —9.2382 
MeV. In Table I the proton pairing energy Eke from the 
SRHBHO and deformed RHB calculations for ?°Mg are 
given. The deformed RHB calculation using the param- 
eter set NL3 with Vj = 380 MeV fm, E&P = 60 MeV 
and the smooth parameter [2 = 5.65 MeV reproduces 
the proton pairing energy from the SRHBHO calculation 
for ?°Mg. These parameters for the pairing are used in 
all the following calculations regardless of whether NL3 


or PK1 is used for the RMF Lagrangian density. 


C. Completeness of the Woods-Saxon basis 


The spherical nucleus ?°?Mg has been investigated as 
the first test of the deformed RHB theory and some re- 
sults were given in Fig. 1 in Ref. [42]. A comparison was 
made between results obtained for ground state prop- 
erties of the spherical nucleus ?°?Mg with the spherical 
RCHB code [29] based on the Runge-Kutta method in 
the radial coordinate r and the new deformed RHB code 
discussed in this manuscript. We summarize these com- 
parison in Table II. In these calculations, the parameter 


~-e NL3, prolate —— NL3, oblate 


42Mg 


FIG. 1. (Color online) Bulk properties of ground state (a) and 
the oblate minimum (b) of “’Mg as functions of the cutoff 
energy E*,. From the lowest to the top panels, the total 
binding energy Eg, the quadrupole deformation 6, and the 
rms radius R are plotted. The parameter set is NL3, the box 


size is Rpox = 20 fm and the step size is Ar = 0.1 fm. 


set NL3, a box of the size Rpox = 4roA!/3 = 13.0 fm 
and a step size Ar = 0.1 fm are used. The surface 6 
pairing force is used with the strength Vo = —374 MeV 
fm° and psat = 0.152 fm. A sharp cutoff is applied on 
the quasiparticle states with E&R = 60 MeV. It is shown 
that when the basis size increases, the total binding en- 
ergy E, the proton pairing energy Ete and the rms 
radius R all converge to the corresponding exact values. 
In practical calculations, Et, may be chosen according 
to the balance between the desired accuracy and the com- 
putational cost. It is concluded [42] that for light nuclei, 
one can safely use E$ = 100 MeV which results in accu- 
racies in the total binding energy and the proton pairing 
energy of about a hundred keV and in the rms radius of 
around 0.002 fm. 

Since we are also interested in drip-line nuclei, next 
we study the dependence of the deformed RHB results 
on the completeness of the Woods-Saxon basis for a very 
neutron rich nucleus. In this subsection we study the 
results with different values of E+,,. For the calculation 
with a Woods-Saxon basis [38], a box of the size Ryox = 
Ar) A!/3 with ro = 1.2 fm is used. In this case Rpox is 
different for different magnesium isotopes, e.g., 13.0 fm 
for 2°Mg and 16.7 fm for *7Mg. In the present work, we 
prefer to use a fixed box size Rpox = 20 fm which is large 
enough for all magnesium isotopes. The mesh size for the 
radial wave function of each Woods-Saxon state is taken 
as 0.1 fm. 

For *#Mg both prolate and oblate minima in the po- 
tential energy surface are searched for and it is found 
that the ground state is prolate. In Fig. 1 the total bind- 
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~-e-- NL3, prolate —=— NL3, oblate 


42 


R (fm) 
R (fm) 


B 
B 


Ep, (MeV) 
Ep (MeV) 


0.1 0.2 0.3 
Ar (fm) 


FIG. 2. (Color online) Bulk properties for *?Mg as in Fig. 1 
but now as functions of the step size Ar. The box size is 
Rox = 20 fm and an energy cutoff is Et, = 100 MeV. 


ing energy Ep, the quadrupole deformation 6, and the 
rms radius R are plotted as functions of Et, for the pro- 
late ground state and for the oblate minimum of 4?Mg, 
respectively. Apparently, when we increase Ef, these 
quantities all converge well. Similar as in the case of 
the spherical nucleus 2°Mg, for light deformed nuclei, the 
cutoff Et, = 100 MeV results in relative accuracies of 
0.5% for the quadrupole deformation, 0.05% for the rms 
radius, and 0.1 % for the total binding energy. 

The box size Rpox = 20 fm and the cutoff energy 
Ex, = 100 MeV are fixed when we investigate the con- 
vergence of the deformed RHB results with respect to 
the mesh size Ar. In Fig. 2 it is shown that when 
the mesh size decreases, the total binding energy Ep, 
the quadrupole deformation 6, and the rms radius R all 
converge well. The difference of the binding energy be- 
tween calculations with Ar = 0.1 fm and Ar = 0.05 
fm is smaller than 0.025 MeV for both minima, which is 
about 0.008% of the total binding energy. When Ar is 
decreased from 0.1 fm to 0.05 fm, the relative changes of 
the quadrupole deformation 8 and the radius R are both 
smaller than 0.01%. 

Figure 3 shows the same quantities as a function of 
the box size Ryox. The relative deviations between the 
rms radius R at Rpox = 20 fm and Rpox = 30 fm are 
about 0.1% for the prolate ground state and 0.01% for 
the oblate minimum. The box size Rpox = 20 fm gives 
also a good accuracy For the quadrupole deformation 8 
and the binding energy. 

In conclusion, in the following calculations, we fix the 
box size at Ryox = 20 fm, the mesh size at Ar = 0.1 
fm, and the cutoff energy for positive energy states in 
the Woods-Saxon basis at Et, = 100 MeV. In each «- 
block, the number of negative energy states in the Dirac 


—— NL3, oblate 


R (fm) 


10 15 20 25 30 10 15 20 25 30 
Rox (fm) Rox (fm) 


FIG. 3. (Color online) Bulk properties for “’Mg as in Fig. 1 
but now as functions of the box size Rpox. The step size is 
Ar = 0.1 fm and the energy cutoff is Et, = 100 MeV. 


sea is the same as that of positive energy states above 
the Dirac gap. The cutoff parameter for A in the expan- 
sion Eq. (14), Amax = 4 and the cutoff parameter for the 
angular quantum number « in the expansion Eq. (17) 
is |kmax| = 10. With these values we do not introduce 
sizable errors. 


IV. RESULTS AND DISCUSSIONS 


In this section, we present results from the deformed 
RHB theory in continuum. We choose magnesium iso- 
topes as examples. After discussing the bulk properties 
of magnesium isotopes, we will focus on the neutron rich 
nucleus 4?Mg. 


A. Bulk properties of magnesium isotopes 


Figure 4 shows the neutron Fermi energy A, and two 
neutron separation energy S2, of magnesium isotopes cal- 
culated with the parameter sets NL3 [81] and PK1 [79]. 
The separation energies are compared with data taken 
from Ref. [84]. Except the different prediction of the 
two-neutron drip line nucleus, the results of the neutron 
Fermi surfaces and two neutron separation energies are 
very similar for both parameter sets. The calculated two 
neutron separation energies Sən of magnesium isotopes 
agree reasonably well with the available experimental val- 
ues except for °7Mg. The large discrepancy in 3?Mg is 
connected to the shape and the shell structure at N = 20 
and will be discussed later. 

Experimentally the nucleus 4°?Mg has been ob- 


1G 20 24 28 32 36 40 44 48 


S>, (MeV) 


16 20 24 28 32 36 40 44 48 
A 


FIG. 4. (Color online) The neutron Fermi energy An (a) and 
two neutron separation energy S2n (b) of magnesium isotopes 
calculated with the parameter sets NL3 and PK1. The data 
of Sən (labeled as “Exp.”) are taken from Ref. [84]. 


served [85]. Theoretically there are several predictions 
on the last bound nucleus in Mg isotopes, e.g., “4Mg 
in the phenomenological finite range droplet model [86], 
40Mg in a macroscopic- microscopic model [87], a RMF 
model with the parameter set NLSH [88], and the Skyrme 
HFB model with the parameter set SLy4 and solved in 
a 3-dimensional Cartesian mesh [52], and +?7Mg from the 
Skyrme HFB model with SLy4 but solved in a trans- 
formed harmonic oscillator basis [37] and the HFB21 
mass table [89]. Therefore the prediction of the two- 
neutron drip line nucleus in Mg isotopes is both model 
and parametrization dependent. In our deformed RHB 
calculations with the parameter set NL3, *°Mg is the 
last nucleus of which the neutron Fermi surface is neg- 
ative and the two neutron separation energy is positive. 
However, with the parameter set PK1, *?Mg is predicted 
to be the last nucleus within the two-neutron drip line. 


The comparison of the quadrupole deformation 8 be- 
tween the theory and the experiment are given in Fig. 5. 
The experimental values of 8 is extracted from the mea- 
sured B(E2 : Of — 27) values and therefore only ab- 
solute values are available [90]. Generally speaking, the 
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FIG. 5. (Color online) The quadrupole deformation param- 
eter 8 of magnesium isotopes calculated with the parameter 
sets NL3 and PK1. The experimental vales (“Exp.”) are 
taken from Ref. [90]. 


ground state quadrupole deformations 8 calculated with 
both parameter sets reproduce the data rather well. Ex- 
ceptions are the nuclei 3*Mg, which turns out to be 
spherical in both models and 3°Mg, which is prolate and 
slightly less deformed than the experiment for PK1 and 
slightly oblate for NL3. In 32Mg, the gap between the 
neutron levels 1d3/2 and 1f7/2 is almost 7 MeV which 
results in a strong closed shell at N = 20. Therefore 
the deformed RHB calculations with both parameter sets 
predict spherical shapes for this nucleus. This also results 
in a large discrepancy from the experiment for the two 
neutron separation energy Sən of 3?Mg as it is seen in 
Fig. 4. Other mean field models predict spherical or al- 
most spherical shapes for 3?°Mg too [52, 88, 91-96]. For 
the isotopes beyond this nucleus with 32 < A < 46 we 
observe large deformations, the so called “island of inver- 
sion” [97-100] which is related to the quenching of the 
N = 20 shell closure. On the mean field level the nucleus 
32Mg does not belong to this island yet. In fact, going 
beyond mean field and calculating the energy surface as 
a function of the deformation parameters one finds that 
this nucleus is a transitional nucleus with an extended 
shoulder reaching to large deformations. This leads in 
GCM calculations with the Gogny force [93] to wavefunc- 
tions with large fluctuations in deformation space and to 
a large B(E2 : OT — 27) value as it is observed in the ex- 
periment [101]. So far it is an open question, why other 
GCM calculations based on Skyrme forces [102] or on 
the relativistic point coupling model PC-F1 [96] cannot 
reproduce this fact. 

Up to *#Mg, the deformed RHB results from the pa- 
rameter set NL3 are very similar to those from the param- 
eter set PK1. Therefore in the following we will mainly 
focus our discussion on the results from PK1. 

In Fig. 6, the root mean square radii for magnesium iso- 
topes are plotted as functions of the neutron number. We 
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FIG. 6. (Color online) The root mean square radii for magne- 
sium isotopes are plotted as functions of the neutron number. 
We display the neutron radius Rn, the proton radius Rp, the 
matter radius R+, and the available data for R+ [51, 103]. The 


roA*/? curve is included to guide the eye. 


display neutron radii Ry, proton radii Rp, matter radii 
R,, the ro A!/? curve with ro = 1 fm, and experimen- 
tal matter radii [51, 103]. The proton radius are almost 
a constant with a very slow increase with increasing N 
due to the neutron-proton coupling included in the mean 
field. With the neutron number increasing, the neutron 
radius R» increases monotonically with an exception at 
32Mg. The neutron radius of 3*Mg is relatively small, 
which is again due to the strong shell effect at N = 20 
in the mean field calculations. It is shown that the de- 
formed RHB results agree well with the experiment for 
the matter radius. The calculated matter radius follows 
roughly the r9.A!/3 curve up to A = 34. From 3ĉMg on, 
the matter radius lies much high above the rọ A‘/* curve. 
This may indicate some exotic structure in these nuclei. 


Figure 7 shows neutron density profiles of even-even 
magnesium isotopes with A > 28 calculated with the pa- 
rameter set PK1. pn,\à=0(r) represents the spherical com- 
ponent of the neutron density distribution (cf. Eq. 14). 
Palri = 0) with ry = y/x? +y? refers to the den- 
sity distribution along the symmetry axis z (@ = 0°) and 
Pn(z = 0,r_) refers to that perpendicular to the symme- 
try axis z (0 = 90°). With increasing A, the spherical 


FIG. 7. (Color online) Neutron density profiles of even-even 
magnesium isotopes with A > 28 calculated with the param- 
eter set PK1. Details are given in the text. 


component of the neutron density distribution py, =0(r) 
changes rapidly at 47Mg. The density distribution along 
the symmetry axis pn(z, r1 = 0) changes abruptly from 
32Mg to 34Mg. This can be understood easily by the 
change in shape in going from the spherical °*Mg to the 
prolate 3*Mg where the density is elongated along the 
z axis. In the direction perpendicular to the symmetry 
axis, the neutron density pn(z = 0,r1) of *7Mg extends 
very far away from the center of the nucleus and a long 
tail emerges, revealing the formation of a halo. 

By comparing pa(z,rı = 0) and pn(z = 0,r_) for 
“Me. it is found that in the tail part, the neutron density 
extends more along the direction perpendicular to the 
symmetry axis. Since this nucleus as a whole is prolate, 
it indicates that the neutron tail has a different shape 
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as the core. This fact is similar to the decoupling of the 
shape of the halo from the shape of the core found for 
tMg in Ref. [43]. Next we will concentrate on 47Mg and 
discuss in details the structure of its ground state. 


B. Ground state of “7Mg 


TABLE III. Properties of “*Mg at the ground state and at 
the oblate minimum derived from deformed RHB calculations 
with the parameter sets NL3 and PK1. The neutron and 
proton Fermi surface Ay and Ap, neutron, proton and total 
quadrupole deformation n, Pp, 6, neutron, proton and to- 
tal radii Rn, Rp, Rt, neutron and proton pairing energies 


EBair, Epai and total binding energy Ep are listed. 
NL3 
An —0.6147 —0.1753} —0.8805  —0.3989 
Xp —24.6731 —23.9050| —24.2695 —22.8118 
Bn —0.3282 0.4155} —0.3299 0.4181 
Bp —0.2426 0.3911] —0.2426 0.3917 
Bs —0.3038 0.4085} —0.3049 0.4105 
Ra 4.0250 4.1077 4.0291 4.0971 
Rp 3.1208 3.1499 3.1393 3.1673 
Ri 3.7888 3.8584 3.7962 3.8544 
EPair| —18.2511 —6.2620| —17.1509 —6.1595 
EBair| —7.0405 0.0000] —6.7639 0.0000 
Eg |—265.4629 —266.4505 | —270.6907 —270.6993 


In the calculations based on the parameter set PK1, 
the chain of Mg isotopes reaches the two-neutron drip 
line at the nucleus #?°Mg. Its properties are summarized 
in Table II. For 4?Mg we find two minima in the en- 
ergy surface as a function of the deformation parameter 
B. The lower one has a prolate shape and corresponds to 
the ground state of 47Mg. The second minimum has an 
oblate shape. From RMF calculations allowing for triax- 
ial deformations [104] we know, however, that the oblate 
minimum is not stable. It forms a saddle point in the 
(3-7) plane and therefore it does not correspond to an 
isomeric state. The ground state is well deformed with a 
quadrupole deformation 6 ~ 0.41, and a very small two 
neutron separation energy So, % 0.22 MeV. The density 
distribution of this weakly bound nucleus has a very long 
tail in the direction perpendicular to the symmetry axis 
(cf. Fig. 7), which indicates the prolate nucleus *7Mg has 
an oblate halo. 

The density distribution in Fig. 8 is decomposed into 
contributions of the oblate “Halo” and of the prolate 
“Core”. Details of this decomposition will be given fur- 
ther down. This indicates the decoupling between the 
deformations of the core and the halo. 

Pairing correlations play a very important role in the 
formation of the halo [27]. For the parameter set PK1 
we find in Table III in the ground state of 4*Mg a van- 
ishing pairing energy for protons and a paring energy 
EPair = —6-26 MeV for the neutrons. For the zero 
range pairing interaction in Eq. (11) only spin singlet 
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FIG. 8. (Color online) Density distributions of the ground 
state of *’Mg with the z axis as the symmetry axis: (a) the 
neutron halo, and (b) the neutron core. 


(S = 0) states and elements diagonal in the quantum 
number p are taken into account in the pairing tensor. 
See appendix F for more details concerning this assump- 
tion. In Fig. 9 we show the components K{*(r) in Eq. 
(E5)and x) (r) in Eq. (E6) of the pairing tensor in the 
ground state of 47Mg for the parameter set PK1. Fig- 
ure 9(b) shows the main component «}*(r) correspond- 
ing to the large components of the Dirac spinor. Com- 
paring Fig. 9(a) and Fig. 9(b) one finds that K\~ (r) is 
smaller by two orders of magnitude than «{*(r). The 
same sign for the quadrupole (A = 2) and the spherical 
(A = 0) components can be understood by the fact that 
the ground state of 4?Mg is prolate in the present calcu- 
lation. The maximum of «{*(r) appears at about 4.8 fm 
indicating that paring in nuclei is a surface effect. The 
hexadecapole components (\ = 4) are much smaller than 
the spherical components (À = 0). 
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FIG. 9. (Color online) The neutron pairing tensor r?«,~ (r) 
(a) and r?x}*(r) (b) with A = 0, 2, and 4 of the ground state 
of Mg from the deformed RHB theory in continuum with 
the parameter set PK1. 


Weakly bound orbitals or those embedded in the con- 
tinuum play a crucial role in the formation of a nuclear 
halo [27, 105, 106]. In order to have an intuitive under- 
standing of the single particle structure, the canonical ba- 
sis is constructed by the method given in Ref. [29]. The 
single particle spectrum around the Fermi level for the 
ground state of 47Mg is shown in Fig. 10. For an axially 
deformed nucleus with spatial reflection symmetry, the 
good quantum numbers of each single particle state in- 
clude the parity 7 and the third component of the angular 
momentum m (labeled by the Nilsson quantum number 
Q in the figures). The occupation probabilities v? in the 
canonical basis have BCS-form [75] and are given by the 
length of the horizontal lines in Fig. 10. To guide the eye 
we also show by a blue dashed line the BCS-formula cal- 
culated with an average gap parameter. The levels close 
to the threshold are labeled by the number i according to 
their energies, and their conserved quantum number Q” 
as well as the main spherical components are given at the 
right hand side. The neutron Fermi level is within the pf 
shell and most of the single particle levels have negative 
parities. Since the chemical potential An ~ —175 keV is 
negative, the corresponding density p(r) is localized and 
the particles occupying the levels in the continuum are 
bound [26]. Since the chemical potential An is close to the 
continuum, orbitals above the threshold have noticeable 
occupations due to the pairing correlations. For instance, 
the occupation probability of the fifth level (Q7 = 3/27) 
is 31.5%. The fourth level Q7 = 1/27 is just below the 
threshold with a single particle energy in the canonical 
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FIG. 10. (Color online) Single neutron levels of ground state 
of ° Mg in the canonical basis as a function of the occupa- 
tion probability v?. The order i, good quantum numbers 97, 
and the main spherical components for orbitals close to the 
threshold are also given. The blue dashed line corresponds to 
the BCS-formula with an average pairing gap. 


basis Ecan = —0.234 MeV and an occupation probability 
of 53.0%. All the other levels below that orbital are well 
bound with £can < —2 MeV. Similar to those of “4Mg in 
Ref. [43], the single neutron levels of 47Mg can be divided 
into two parts, the deeply bound levels (€c¢an < —2 MeV) 
corresponding to the “core”, and the remaining weakly 
bound levels close to the threshold (€can > —0.3 MeV) 
and in the continuum corresponding to the “halo”. 

We have already seen in Fig. 8 that the core is pro- 
late and the halo is oblate. According to Eq. (14) the 
density distributions of the core and of the halo are de- 
composed into spherical (A = 0), quadrupole (A = 2), 
and hexadecapole (A = 4) components in Fig. 11. The 
quadrupole component of the core turns out to be posi- 
tive, which is consistent with the prolate shape of 4?7Mg in 
the ground state. However, for the halo, the quadrupole 
component is mainly negative, which means the halo has 
an oblate shape. This explains the decoupling between 
the quadrupole deformations of the core and the halo. 
We also find in Fig. 11 that the spherical component is 
absolutely the main part of the density distribution for 
both the core and the halo, and that the hexadecapole 
component in the density distribution of the neutron halo 
is also noticeable. 

In order to study the formation mechanism of the halo 
in more detail, we show in Fig. 12(a) the main (spher- 
ical) components ph, y=o Of the density distribution for 
the weakly bound neutron orbitals i. Figure 12(b) gives 
the ratio of these spherical components pi y=9 to the 
spherical component of the total neutron density pn,,=0.- 
One can clearly see that far away from the center, the 
main contribution comes from the 4th and 5th levels. 
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FIG. 11. (Color online) Decomposition of the neutron density 
of the ground state of #*Mg into spherical (A = 0), quadrupole 
(A = 2), and hexadecapole (A = 4) components for the halo 
(a) and the core (b). 


Almost 80% of the total density distribution in the tail 
part comes from these two levels which are close to the 
Fermi surface. Level 7 is embedded in the continuum and 
gives also some contribution to the tail of the total den- 
sity distribution. However, the occupation probability of 
this level is just 5.7%, so its contribution is very small. 
The occupation probability of level 6 is 7.9%, a bit lager 
than that of level 7. But there is almost no contribution 
to the tail of total density from this level. By examin- 
ing the spherical Woods-Saxon components, it is found 
that the main component of level 6 is 1f7/2. The large 
centrifugal barrier of f states with l = 3 hinders strongly 
its spatial extension. For level 7, about 31.3% contribu- 
tion comes from 2pı/2 with a small centrifugal barrier 
and therefore the density can extend far away from the 
center of the nucleus. 


As it is shown in Fig. 12, the halo is mainly formed by 
level 4 and level 5 with occupation probabilities of 53.0% 
and 31.5% respectively. Having in mind the degeneracy 
2 for each single particle level, the occupation number of 
these two orbitals is about 1.7. If we decompose the de- 
formed wave functions of these two orbitals in the spher- 
ical Woods-Saxon basis, it turns out that in both cases 
the major part comes from p waves, as indicated on the 
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FIG. 12. (Color online) Spherical components of neutron den- 
sity distributions of ground state of “Mg: (a) the total den- 
sity and its decomposition into the core and the halo and 
contributions from several neutron orbitals around the Fermi 
level; (b) relative contributions of these neutron orbitals to 
the total neutron density. 


right-hand side of Fig. 10. For level 4 (Q7 = 1/27), the 
probability of 2p3/2, 1fs/2, and 2p;/2 are 37.0%, 32.3%, 
and 21.2% respectively. For level 5 (Q7 = 3/27), 2p372 
is the dominant component with a probability of 78.6%. 
The low centrifugal barrier for p waves gives rise to the 
formation of the halo. 

The shape of the halo originates from the intrin- 
sic structure of the weakly bound or continuum or- 
bitals [43, 45]. As discussed before, for the ground 
state of “?Mg, the halo is mainly formed by level 4 
and level 5. We know that the angular distribution of 
|¥10(0, d)|? x cos? 0 with a projection of the orbital angu- 
lar momentum on the symmetry axis A = 0 is prolate and 
that of |¥i41(0,¢)|? œ sin? @ with A = 1 is oblate [45]. 
For level 4 (Q7 = 1/27), A could be 0 or 1 since the 
third component of total spin is 1/2. However, it turns 
out that the A = 0 component dominates which results in 
an oblate shape. For level 5, since the third component of 
the total spin is 3/2, A can only be 1, which corresponds 
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to an oblate shape too. Therefore in 4?Mg the shape of 
the halo is oblate and decouples from the prolate core. 


V. SUMMARY 


A deformed relativistic Hartree Bogoliubov theory in 
continuum is developed in order to describe deformation 
effects in exotic nuclei allowing for halo structures. The 
deformed RHB equations are solved in a Woods-Saxon 
basis where the radial wave functions have a proper 
asymptotic behavior at large distance from the nuclear 
center. This is crucial for the formation of a halo. The 
formalism and the numerical details of the deformed 
RHB theory are presented. Routine checks are made 
including convergence studies of the deformed RHB re- 
sults concerning the mesh size, the box size and the size 
of the Woods-Saxon basis. The results are compared for 
spherical nuclei with solutions of the 1D continuum RHB 
equations in the radial coordinate r based on the Runge- 
Kutta method. 

The deformed RHB theory in continuum is applied to 
study the chain of magnesium isotopes with the param- 
eter sets NL3 and PK1 of the Lagrangian. Except for 
the different prediction of the two-neutron drip line nu- 
cleus, the results of neutron Fermi surfaces and two neu- 
tron separation energies are very similar for both parame- 
ter sets. The calculated two neutron separation energies 
So, of magnesium isotopes agree reasonably well with 
the available experimental values except for 9*Mg, a well 
known problem connected with the shape and the shell 
structure at N = 20. For 3?Mg, the gap between the neu- 
tron levels 1d3/2 and 1f7/2 is almost 7 MeV which results 
in a strong shell closure at N = 20. The nuclear radii are 
also investigated, the deformed RHB results agree well 
with the experiment for matter radii. The proton radius 
is almost a constant with a very slow increase with in- 
creasing N due to the neutron-proton coupling included 
in the mean field. A sharp increase in the neutron radius 
is observed at *?Mg. 

Detailed results are shown for the two-neutron drip line 
nucleus 4?Mg with the parameter set PK1, which is well 
deformed. The ground state of 4?Mg is prolate, however, 
it has an oblate neutron halo. By examining in detail the 
density distributions, the pairing tensor, and the single 


[1] J. W. Xia, W. L. Zhan, B. W. Wei, Y. J. Yuan, M. T. 
Song, W. Z. Zhang, X. D. Yang, P. Yuan, D. Q. Gao, 
H. W. Zhao, X. T. Yang, G. Q. Xiao, K. T. Man, J. R. 
Dang, X. H. Cai, Y. F. Wang, J. Y. Tang, W. M. Qiao, 
Y. N. Rao, Y. He, L. Z. Mao, and Z. Z. Zhou, Nucl. 
Instrum. Methods Phys. Res. A 488, 11 (2002). 

[2] W. Zhan, H. Xu, G. Xiao, J. Xia, H. Zhao, and Y. Yuan, 
Nucl. Phys. A 834, 694c (2010). 

[3] C. Sturm, B. Sharkov, and H. Stöcker, Nucl. Phys. A 
834, 682c (2010). 


12 


particle levels in the canonical basis in the deformed nu- 
cleus 4?Mg, it can be understood, why the shape of the 
neutron halo decouples from that of the core. It is shown 
that the existence and the deformation of a possible neu- 
tron halo depends essentially on the quantum numbers 
of the main components of the single-particle orbits in 
the vicinity of the Fermi surface and the shape of their 
single-particle density distributions. 

In stable nuclei, there are situations that the levels of 
valence nucleons are sometimes also well separated from 
the core. It is, however, a difficult question, whether 
there exists cases of such a decoupling of shapes as we 
have seen here in the case of loosely bound valence orbits 
close to the continuum limit, because in stable nuclei 
even the valence nucleons are well bound in the average 
potential. 

We can conclude that spherical and deformed relativis- 
tic Hartree Bogoliubov theory in continuum is a very 
powerful tool providing a proper description of exotic 
nuclei including halo phenomena, because it takes into 
account in a self-consistent and microscopic way polar- 
ization effects, shape changes of individual orbitals, pair- 
ing correlations and the coupling to the continuum with 
proper boundary conditions. 


ACKNOWLEDGMENTS 


This work has been supported in part by the Natu- 
ral Science Foundation of China (10875157, 10975100, 
10979066, 11105005, 11175002, and 11175252), by the 
Major State Basic Research Development Program of 
China (2007CB815000), by the Knowledge Innovation 
Project of Chinese Academy of Sciences (KJCX2-EW- 
N01 and KJCX2-YW-N32), and by the DFG cluster 
of excellence “Origin and Structure of the Universe” 
(www.universe-cluster.de). Part of the results described 
in this paper is obtained on the ScGrid of Supercom- 
puting Center, Computer Network Information Center of 
Chinese Academy of Sciences. One of the authors (P.R.) 
would like to express his gratitude to J. Meng for the 
kind hospitality extended to him at the Peking Univer- 
sity. Helpful discussions with N. V. Giai, B. N. Lu, Z. Y. 
Ma, N. Sandulescu, J. Terasaki, D. Vretenar, S. J. Wang, 
and S. Yamaji are gratefully acknowledged. 


4] S. Gales, Nucl. Phys. A 834, 717c (2010). 


5| T. Motobayashi and Y. Yano, Nucl. Phys. News 17, 
5 (2007); T. Motobayashi, Nucl. Phys. A 834, 707c 
(2010). 


6] M. Thoennessen, Nucl. Phys. A 834, 688c (2010). 

S. Choi, “KoRIA project - RI accelerator in Korea,” In- 
vited talk given in International Symposium on Nuclear 
Physics in Asia, 14-15 October, 2010, Beihang Univer- 
sity, Beijing. 

8] C. A. Bertulani, M. S. Hussein, 


~ 


and G. Münzenberg, 


201708.00252v1 


chinaXiv 


Physics of Radioactive Beams (Nova Science Publishers, 
Inc., 2001). 

[9] A. C. Miiller and B. M. Sherrill, Annu. Rev. Nucl. Part. 

Sci. 43, 529 (1993). 

10] I. Tanihata, Prog. Part. Nucl. Phys. 35, 505 (1995). 

11] P. G. Hansen, A. S. Jensen, and B. Jonson, Annu. Rev. 

Nucl. Part. Sci. 45, 591 (1995). 

12] R. F. Casten and B. M. Sherrill, Prog. Part. Nucl. Phys. 

45, S171 (2000). 

13] B. Johnson, Phys. Rep. 389, 1 (2004). 

14] A. S. Jensen, K. Riisager, D. V. Fedorov, and E. Gar- 

rido, Rev. Mod. Phys. 76, 215 (2004). 

15] S. N. Ershov, L. V. Grigorenko, J. S. Vaagen, and M. V. 

Zhukov, J. Phys. G: Nucl. Phys. 37, 064026 (2010). 

16] Z.-X. Cao and Y.-L. Ye, Sci. China-Phys. Mech. Astron. 

54 (Suppl. 1), s1 (2011). 

17| I. Tanihata, H. Hamagaki, O. Hashimoto, Y. Shida, 
N. Yoshikawa, K. Sugimoto, O. Yamakawa, 
T. Kobayashi, and N. Takahashi, Phys. Rev. Lett. 55, 
2676 (1985). 

[18] T. Minamisono, T. Ohtsubo, I. Minami, S. Fukuda, 
A. Kitagawa, M. Fukuda, K. Matsuta, Y. Nojiri, 
S. Takeda, H. Sagawa, and H. Kitagawa, Phys. Rev. 
Lett. 69, 2058 (1992). 

[19] W. Schwab, H. Geissel, H. Lenske, K. H. Behr, 
A. Brünle, K. Burkard, H. Irnich, T. Kobayashi, 
G. Kraus, A. Magel, G. Münzenberg, F. Nickel, K. Riis- 
ager, C. Scheidenberger, B. M. Sherrill, T. Suzuki, and 
B. Voss, Z. Phys. A 350, 283 (1995). 

[20] A. Ozawa, T. Kobayashi, T. Suzuki, K. Yoshida, and 
I. Tanihata, Phys. Rev. Lett. 84, 5493 (2000). 

[21] P. Adrich, A. Klimkiewicz, M. Fallot, K. Boretzky, 
T. Aumann, D. Cortina-Gil, U. D. Pramanik, T. W. 
Elze, H. Emling, H. Geissel, M. Hellstrom, K. L. Jones, 
J. V. Kratz, R. Kulessa, Y. Leifels, C. Nociforo, R. Palit, 
H. Simon, G. Surowka, K. Summerer, and W. Walus, 
Phys. Rev. Lett. 95, 132501 (2005). 

22] L.-g. Cao and Z.-y. Ma, Phys. Rev. C 66, 024311 (2002). 

23] J. Dobaczewski, N. Michel, W. Nazarewicz, M. Plosza- 

jezak, and J. Rotureau, Prog. Part. Nucl. Phys. 59, 432 

(2007). 

24] J. C. Pei, A. T. Kruppa, and W. Nazarewicz, Phys. 

Rev. C 84, 024311 (2011). 

25] A. Bulgac, “Hartree-Fock-Bogoliubov approximation 

for finite systems,” IPNE FT-194-1980, Bucharest 

(arXiv: nucl-th/9907088) (1980). 

26] J. Dobaczewski, H. Flocard, and J. Treiner, Nucl. Phys. 
A 422, 103 (1984); J. Dobaczewski, W. Nazarewicz, 
T. R. Werner, J. F. Berger, C. R. Chinn, and 
J. Dechargé, Phys. Rev. C 53, 2809 (1996). 

[27] J. Meng and P. Ring, Phys. Rev. Lett. 77, 3963 (1996). 

[28] W. Pöschl, D. Vretenar, G. A. Lalazissis, and P. Ring, 
Phys. Rev. Lett. 79, 3841 (1997); G. Lalazissis, 
D. Vretenar, W. Pöschl, and P. Ring, Phys. Lett. B 
418, 7 (1998). 

[29] J. Meng, Nucl. Phys. A 635, 3 (1998). 

[30] N. Sandulescu, N. Van Giai, and R. J. Liotta, Phys. 
Rev. C 61, 061301(R) (2000); N. Sandulescu, L. S. 
Geng, H. Toki, and G. C. Hillhouse, Phys. Rev. C 68, 
054323 (2003). 

[31] S.-S. Zhang, E.-G. Zhao, and S.-G. Zhou, “Theoretical 
study of the two-proton halo candidate "Ne including 
contributions from a resonant continuum and pairing 
correlations,” arXiv: 1105.0504 [nucl-th] (2011), sub- 


13 


mitted to Nucl. Phys. A. 

[32] W. H. Press, S. A. Teukolsky, W. T. Vetterling, and 
B. P. Flannery, Numerical Recipes in Fortran 77: the 
Art of Scientific Computing (Vol. 1 of Fortran Nu- 
merical Recipes), 2nd ed. (Cambridge University Press, 
1992). 

[33] W. Pöschl, D. Vretenar, A. Rummel, and P. Ring, 
Comput. Phys. Commun. 101, 75 (1997); W. Péschl, 
D. Vretenar, and P. Ring, Comput. Phys. Commun. 
103, 217 (1997). 

34] D. Vautherin, Phys. Rev. C 7, 296 (1973). 

35] D. Gogny, Nucl. Phys. A 237, 399 (1975); J. Decharge 
and D. Gogny, Phys. Rev. C 21, 1568 (1980). 

36] Y. K. Gambhir, P. Ring, and A. Thimet, Ann. Phys. 
198, 132 (1990); Y. K. Gambhir and P. Ring, Mod. 
Phys. Lett. A 8, 787 (1993). 

37] M. V. Stoitsov, J. Dobaczewski, P. Ring, and S. Pit- 
tel, Phys. Rev. C 61, 034311 (2000); M. V. Stoitsov, 
J. Dobaczewski, W. Nazarewicz, S. Pittel, and D. J. 
Dean, Phys. Rev. C 68, 054312 (2003). 

38] S.-G. Zhou, J. Meng, and P. Ring, Phys. Rev. C 68, 

034323 (2003). 

39] N. Schunck and J. L. Egido, Phys. Rev. C 77, 011301 (R) 

(2008); Phys. Rev. C 78, 064305 (2008). 

40] W. H. Long, P. Ring, N. V. Giai, and J. Meng, Phys. 

Rev. C 81, 024308 (2010). 

41] S.-G. Zhou, J. Meng, and P. Ring, in Nuclear Physics 

Trends, AIP Conf. Proc., Vol. 865, edited by Y.-G. Ma 

and A. Ozawa (AIP, 2006) pp. 90-95. 

42] S.-G. Zhou, J. Meng, and P. Ring, in Physics of Unsta- 
ble Nuclei, edited by D. T. Khoa, P. Egelhof, S. Gales, 
N. Van Giai, and T. Motobayashi (World Scientific, 
2008) pp. 402-408, proceedings of the International 
Symposium on Physics of Unstable Nuclei, July 3-7, 
2007, Hoi An, Vietnam, arXiv: 0803.1376v1 [nucl-th]. 

[43] S.-G. Zhou, J. Meng, P. Ring, and E.-G. Zhao, 
Phys. Rev. C 82, 011301(R) (2010), arXiv:0909.1600v3 
[nucl-th]; J. Phys: Conf. Ser. 312, 092067 (2011), 
arXiv:1101.3158v1 [nucl-th]. 

44] X. Li and P.-H. Heenen, Phys. Rev. C 54, 1617 (1996). 

45| T. Misu, W. Nazarewicz, and S. Åberg, Nucl. Phys. A 

614, 44 (1997). 

46] L. Guo, E.-G. Zhao, and F. Sakata, Commun. Theor. 

Phys. 40, 573 (2003). 

47| J. Meng, H. F. Lü, S. Q. Zhang, and S. G. Zhou, Nucl. 

Phys. A 722, 366c (2003). 

48] F. M. Nunes, Nucl. Phys. A 757, 349 (2005). 

49| J. C. Pei, F. R. Xu, and P. D. Stevenson, Nucl. Phys. 

A 765, 29 (2006). 

50] T. Nakamura, N. Kobayashi, Y. Kondo, Y. Satou, 
N. Aoi, H. Baba, S. Deguchi, N. Fukuda, J. Gibelin, 
N. Inabe, M. Ishihara, D. Kameda, Y. Kawada, 
T. Kubo, K. Kusaka, A. Mengoni, T. Motobayashi, 
T. Ohnishi, M. Ohtake, N. A. Orr, H. Otsu, T. Ot- 
suka, A. Saito, H. Sakurai, S. Shimoura, T. Sumikama, 
and H. Takeda, Phys. Rev. Lett. 103, 262501 (2009). 

[51] R. Kanungo, A. Prochazka, W. Horiuchi, C. Nociforo, 
T. Aumann, D. Boutin, D. Cortina-Gil, B. Davids, 
M. Diakaki, F. Farinon, H. Geissel, R. Gernhaeuser, 
J. Gerl, R. Janik, B. Jonson, B. Kindler, R. Knoebel, 
R. Kruecken, M. Lantz, H. Lenske, Y. Litvinov, B. Lom- 
mel, K. Mahata, P. Maierbeck, A. Musumarra, T. Nils- 
son, C. Perro, C. Scheidenberger, B. Sitar, P. Strmen, 
B. Sun, Y. Suzuki, I. Szarka, I. Tanihata, Y. Utsuno, 


201708.00252v1 


chinaXiv 


[52] 


53 
54 
55 
56 
57 
58 
59 


60 


61 


62 


63 


64 


65 


66 


67 
68 


69 
70 


71 


72 


73 


74 
75 


76 


T 
78 


79 


80 


81 


82 
83 


H. Weick, and M. Winkler, Phys. Rev. C 83, 021302(R) 
(2011). 

J. Terasaki, P. H. Heenen, H. Flocard, and P. Bonche, 
Nucl. Phys. A 600, 371 (1996); J. Terasaki, H. Flocard, 
P. H. Heenen, and P. Bonche, Nucl. Phys. A 621, 706 
(1997). 

M. Yamagami, K. Matsuyanagi, and M. Matsuo, Nucl. 
Phys. A 693, 579 (2001). 

E. Teran, V. E. Oberacker, and A. S. Umar, Phys. Rev. 
C 67, 064314 (2003). 

V. E. Oberacker, A. S. Umar, E. Terán, 
A. Blazkiewicz, Phys. Rev. C 68, 064302 (2003). 
J. C. Pei, M. V. Stoitsov, G. I. Fann, W. Nazarewicz, 
N. Schunck, and F. R. Xu, Phys. Rev. C 78, 064306 
(2008). 

N. Tajima, Phys. Rev. C 69, 034305 (2004). 

P.-G. Reinhard, M. Bender, K. Rutz, and J. Maruhn, 
Z. Phys. A 358, 277 (1997). 

H. Nakada, Nucl. Phys. A 808, 47 (2008). 

H. Oba and M. Matsuo, Phys. Rev. C 80, 024301 (2009); 
Y. Zhang, M. Matsuo, and J. Meng, Phys. Rev. C 83, 
054301 (2011). 

D. Vretenar, G. A. Lalazissis, and P. Ring, Phys. Rev. 
Lett. 82, 4595 (1999). 

G. A. Lalazissis, D. Vretenar, P. Ring, M. Stoitsov, and 
L. M. Robledo, Phys. Rev. C 60, 014310 (1999); G. A. 
Lalazissis, D. Vretenar, and P. Ring, Nucl. Phys. A 
650, 133 (1999). 

T. Niksic, D. Vretenar, P. Ring, and G. A. Lalazissis, 
Phys. Rev. C 65, 054320 (2002). 

T. Niksic, D. Vretenar, G. A. Lalazissis, 
Phys. Rev. C 69, 047301 (2004). 

J.-K. Zhang and D. S. Onley, Phys. Lett. B 209, 145 
(1988); Nucl. Phys. A 526, 245 (1991). 

H.-P. Duerr and E. Teller, Phys. Rev. 101, 494 (1956); 
H.-P. Duerr, Phys. Rev. 103, 469 (1956). 

J. D. Walecka, Ann. Phys. 83, 491 (1974). 

B. D. Serot and J. D. Walecka, Adv. Nucl. Phys. 16, 1 
(1986). 

P. G. Reinhard, Rep. Prog. Phys. 52, 439 (1989). 

P. Ring, Prog. Part. Nucl. Phys. 37, 193 (1996); Prog. 
Part. Nucl. Phys. 46, 165 (2001). 

D. Vretenar, A. Afanasjev, G. Lalazissis, and P. Ring, 
Phys. Rep. 409, 101 (2005). 

J. Meng, H. Toki, S. G. Zhou, S. Q. Zhang, W. H. Long, 
and L. S. Geng, Prog. Part. Nucl. Phys. 57, 470 (2006). 
J. Boguta and A. R. Bodmer, Nucl. Phys. A 292, 413 
(1977). 

H. Kucharek and P. Ring, Z. Phys. A 339, 23 (1991). 
P. Ring and P. Schuck, The Nuclear Many-Body Prob- 
lem (Springer, 1980). 

C. E. Price and G. E. Walker, Phys. Rev. C 36, 354 
(1987). 

W. Koepf and P. Ring, Z. Phys. A 339, 81 (1991). 

P. Bonche, H. Flocard, P. H. Heenen, S. J. Krieger, and 
M. S. Weiss, Nucl. Phys. A 443, 39 (1985). 

W. Long, J. Meng, N. V. Giai, and S.-G. Zhou, Phys. 
Rev. C 69, 034319 (2004). 

P.-W. Zhao, B.-Y. Sun, and J. Meng, Chin. Phys. Lett. 
26, 112102 (2009). 

G. A. Lalazissis, J. Konig, and P. Ring, Phys. Rev. C 
55, 540 (1997). 

M. Serra and P. Ring, Phys. Rev. C 65, 064324 (2002). 
J. F. Berger, M. Girod, and D. Gogny, Nucl. Phys. A 


and 


and P. Ring, 


[84] 


[85] 


86 


87 
88 


89 
90 
91 
92 
93 
94 


95 


96 


[97] 


[98] 


[99] 


[100] 


[101] 


14 


428, 23 (1984). 

G. Audi, A. H. Wapstra, and C. Thibault, Nucl. Phys. 
A 729, 337 (2003). 

T. Baumann, A. M. Amthor, D. Bazin, B. A. Brown, 
C. M. F. III, A. Gade, T. N. Ginter, M. Hausmann, 
M. Matos, D. J. Morrissey, M. Portillo, A. Schiller, 
B. M. Sherrill, A. Stolz, O. B. Tarasov, and M. Thoen- 
nessen, Nature 449, 1022 (2007). 

P. Möller, J. R. Nix, W. D. Myers, and W. J. Swiatecki, 
At. Data Nucl. Data Tables 59, 185 (1995). 

Q. Zhi and Z. Ren, Phys. Lett. B 638, 166 (2006). 

Z. Ren, Z. Y. Zhu, Y. H. Cai, and G. Xu, Phys. Lett. 
B 380, 241 (1996). 

S. Goriely, N. Chamel, and J. M. Pearson, Phys. Rev. 
C 82, 035804 (2010), and the HFB21 mass table. 

S. Raman, C. W. Nestor, and P. Tikkanen, At. Data 
Nucl. Data Tables 78, 1 (2001). 

S. K. Patra and C. R. Praharaj, Phys. Lett. B 273, 13 
(1991). 

G. Lalazissis, S. Raman, and P. Ring, At. Data Nucl. 
Data Tables 71, 1 (1999). 

R. R. Rodriguez-Guzman, J. L. Egido, and L. M. Rob- 
ledo, Phys. Rev. C 62, 054319 (2000). 

R. Rodriguez-Guzman, J. Egido, and L. M. Robledo, 
Nucl. Phys. A 709, 201 (2002). 

T. Niksic, D. Vretenar, and P. Ring, Phys. Rev. C 73, 
034308 (2006). 

J. M. Yao, J. Meng, P. Ring, and D. Vretenar, Phys. 
Rev. C 81, 044311 (2010); J. M. Yao, H. Mei, H. Chen, 
J. Meng, P. Ring, and D. Vretenar, Phys. Rev. C 83, 
014308 (2011). 

O. Sorlin and M.-G. Porquet, Prog. Part. Nucl. Phys. 
61, 602 (2008). 

V. Tripathi, S. L. Tabor, P. F. Mantica, Y. Utsuno, 
P. Bender, J. Cook, C. R. Hoffman, S. Lee, T. Otsuka, 
J. Pereira, M. Perry, K. Pepper, J. S. Pinter, J. Stoker, 
A. Volya, and D. Weisshaar, Phys. Rev. Lett. 101, 
142504 (2008). 

P. Doornenbal, H. Scheit, N. Aoi, S. Takeuchi, K. Li, 
E. Takeshita, H. Wang, H. Baba, S. Deguchi, N. Fukuda, 
H. Geissel, R. Gernhauser, J. Gibelin, I. Hachiuma, 
Y. Hara, C. Hinke, N. Inabe, K. Itahashi, S. Itoh, 
D. Kameda, S. Kanno, Y. Kawada, N. Kobayashi, 
Y. Kondo, R. Krucken, T. Kubo, T. Kuboki, K. Kusaka, 
M. Lantz, S. Michimasa, T. Motobayashi, T. Naka- 
mura, T. Nakao, K. Namihira, S. Nishimura, T. Ohnishi, 
M. Ohtake, N. A. Orr, H. Otsu, K. Ozeki, Y. Satou, 
S. Shimoura, T. Sumikama, M. Takechi, H. Takeda, 
K. N. Tanaka, K. Tanaka, Y. Togano, M. Winkler, 
Y. Yanagisawa, K. Yoneda, A. Yoshida, K. Yoshida, 
and H. Sakurai, Phys. Rev. Lett. 103, 032501 (2009). 
K. Wimmer, T. Kroell, R. Kruecken, V. Bildstein, 
R. Gernhaeuser, B. Bastin, N. Bree, J. Diriken, 
P. Van Duppen, M. Huyse, N. Patronis, P. Vermae- 
len, D. Voulot, J. Van de Walle, F. Wenander, L. M. 
Fraile, R. Chapman, B. Hadinia, R. Orlandi, J. F. 
Smith, R. Lutter, P. G. Thirolf, M. Labiche, A. Blazhev, 
M. Kalkuehler, P. Reiter, M. Seidlitz, N. Warr, A. O. 
Macchiavelli, H. B. Jeppesen, E. Fiori, G. Georgiev, 
G. Schrieder, S. Das Gupta, G. Lo Bianco, S. Nardelli, 
J. Butterworth, J. Johansen, and K. Riisager, Phys. 
Rev. Lett. 105, 252501 (2010). 

J. A. Church, C. M. Campbell, D.-C. Dinca, J. En- 
ders, A. Gade, T. Glasmacher, Z. Hu, R. V. F. Janssens, 


201708.00252v1 


chinaXiv 


W. F. Mueller, H. Olliver, B. C. Perry, L. A. Riley, and 
K. L. Yurkewicz, Phys. Rev. C 72, 054320 (2005). 
[102] P.-H. Heenen, P. Bonche, S. Cwiok, W. Nazarewicz, and 
A. Valor, RIKEN Rev. 26, 31 (2000). 
[103] T. Suzuki, H. Geissel, O. Bochkarev, L. Chulkov, 
M. Golovkov, N. Fukunishi, D. Hirata, H. Ir- 
nich, Z. Janas, H. Keller, T. Kobayashi, G. Kraus, 
G. Müzenberg, S. Neumaier, F. Nickel, A. Ozawa, 
A. Piechaczeck, E. Roeckl, W. Schwab, K. Sümerer, 
K. Yoshida, and I. Tanihata, Nucl. Phys. A 630, 661 
(1998). 
B.-N. Lu, E.-G. Zhao, and S.-G. Zhou, Phys. Rev. C 84, 
014328 (2011), arXiv:1104.4638v1 [nucl-th]; Phys. Rev. 
C 85, 011301 (R) (2012), arXiv: 1110.6769v2 [nucl-th]. 
J. Meng, W. Pöschl, and P. Ring, Z. Phys. A 358, 123 
(1997). 
J. Meng and P. Ring, Phys. Rev. Lett. 80, 460 (1998). 
A. R. Edmonds, Angular Momentum in Quantum Me- 
chanics (University Press, Princeton, 1957). 
R. E. Peierls and J. Yoccoz, Proc. Phys. Soc. A 70, 381 
(1957). 
M. Bender, K. Rutz, P.-G. Reinhard, and J. A. Maruhn, 
Eur. Phys. J. A 7, 467 (2000). 


104 


105 


106 
107 


108 


109 


Appendix A: Spherical spinors in coordinate space 


In this work we use three different representations of 
the wave functions. The starting point is the coordinate 
space representation x = (rsp), where s is the spin coor- 
dinate and p describes large (or upper) (p = 1 or p = +) 
and small (or lower) (p = 2 or p = —) components. The 
second basis is a discrete basis of spherical Dirac spinors 
jnm) which is obtained by the diagonalization of the 
spherical Dirac Hamiltonian with fields of Woods-Saxon 
shape. This basis is called the Woods-Saxon basis in the 
following. In this basis the RHB equations is solved and 
the solutions form a basis of quasi-particle states labeled 
by |&). The Dirac spinors of the Woods-Saxon basis are 
represented in coordinate space as 


(rsp|nkMm) = bnam(Tsp) = gp Rol) 


Yum (9,8), 

(A1) 
The orbital angular momenta of these components are 
l(p = 1) = j + $sign(«) and I(p = 2) = j — įsign(«). 
Ran) = Gnelr), Ras(r,2) = Fne(r) are the radial 
wave functions, and Y,.,, are the spinor spherical har- 
monics 


1 ; 
Yim) > Olgmslmi|jm)Yim (Q)X 4m, (A2) 
M1, Ms 
The time reversal state reads 


Dace (rsp) = (—1)PHP)+5-™ 4g, om (rsp). (A3) 


ee ee oe aa 
2 —A*(rsp,r's'p!) —hp(rsp,r's'p!) +) ve(risip!)) T N Velrsp) 
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These basis functions are obtained from the solution 
of a Dirac equation with spherical potentials of Woods- 
Saxon-shape [77] 


hD =a-p+B |v +50 (| +VO(r), (A4) 
on a mesh in r-space using the Runge-Kutta method. 
For each «x we have eigenstates with positive and nega- 
tive eigenvalues €,,,, and for completeness of the basis the 
sum over nw has to include states with positive eigenval- 
ues and those with negative eigenvalues [38]. This has 
nothing to do with the no-sea approximation which is 
applied in the final quasiparticle basis where the sums 
over k in Eq. (9) runs only over solutions with positive 
single particle energies. 

Since the RHB equation (19) has to be solved in this 
basis one has to evaluate matrix elements of the form 


(nkm|hp|n'«’m) and (nam|A|n'k'’m). (A5) 
In order to simplify the calculations, the integrations over 
the angles are carried out analytically using well known 
angular momentum coupling techniques and only the ra- 
dial integrals are calculated numerically. For local po- 
tentials we need the following products of basis wave- 
functions 


5 Onrm. (TSP) Oh etm (rsp). (A6) 


Following Eq. (14) they are expanded in terms of Legen- 
dre polynomials. For the coefficient of rank ÀA depending 
only on the radius r we find 


Z Rnr(r, p) Raie (r, p) (km|Pi|k'm) 
x T r i 


5 Onini Oniki 
l (A7) 


The angular matrix elements (sm|P,|k’'m) can be de- 
rived with the help of the Wigner-Eckart theorem [107]. 
For even values of l + A +l we find 


. ag ag 
impe = oaa (FF) (404), 
(A8) 
where j = 27 +1. For odd values of l + A +l these 
matrix elements vanish. 


Appendix B: Matrix elements of the DRHB 
Hamiltonian 


The Dirac Hartree-Bogoliubov equations [74] read in 
coordinate space 


r sp 
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where Ey is the quasiparticle energy and A the chemical 
potential. On the Hartree level the Dirac Hamiltonian is 
local 
hp(rsp,r’s'p') = hp(r, sp, s'p')0(r — r’). (B2) 
For the zero range pairing force in Eq. (11) which projects 
onto the S = 0 part of the pairing density, the pairing 
field is local too and does not depend on the spin variables 
A(rp,r'p!) = Spy Alrp)o(r =r"). (B3) 
In this work we restrict ourselves on pairing fields diag- 
onal in the quantum number p (see Appendix F). These 
equations of motions are solved by expanding the spinors 
U; and V;, in terms of a Woods-Saxon basis of Dirac 


16 


spinors Ynem(rs) in Eq. (18) with positive and negative 
single particle energies €,,;. 

For the self-consistent solution of the Dirac equation 
(5) with deformed potentials of axial symmetry, we ex- 
pand the potentials S(r) and V(r) in terms of the Leg- 
endre polynomials as in Eq. (14). The deformed Dirac 
Hamiltonian hp is divided into two parts, the spherical 
Woods-Saxon Hamiltonian hO of Eq. (A4) and the de- 
formed rest 

hp =t ES) 


) + V(r] PA(Q), (BA) 


with S$ = So — 5, W =V, —V©), and S$ = Sy, and 
VX = Vy for à > 0. Using Eq. (A7) the matrix elements 
of the Dirac Hamiltonian read, 


(nk|Ap|n's’) = en Onn’ On! + SD trm|Palatm) f dr [Gnr (V(r) + SY (7)) Gin (7) 


A 


The integral in the pairing matrix element 
(nkm|A|n'k'm) contains the time reversal basis function. 
Since the pairing interaction Eq. (11) projects onto the 
S = 0 we have to couple the product @,,4,(178)@n'x'm(T8) 
to spin S = 0 and find 


p Prkl S 


Pn’ pintie a): 


(B6) 


z= Pnrmls Jn’ K! mf 


`z 
1 
NI 


Using again Eq. (A7) one finds 


(nk Attn’ K’) = S iem|Pale’m) f drGus Aa) Gnw 
X 
(B7) 
and 


(nk ATT nk’) = Siem Paje’) | drFas Aat) Fw 
A 
(B8) 
where the potentials A, (r) will be given in Appendix E. 


Fnr lr) (V(r) — S\(r)) Fai 
————_ 


Appendix C: Calculation of the densities 


In order to determine the self-consistent fields in the 
next step of the iteration we first have to determine 
the densities. Starting from the expansion coefficients 
un T and ae obtained through the diagonalization 
of the RHB matrix (19) we find the density matrix in the 
Woods-Saxon basis 


(m) = (m)* (m) 
Pini) (n'k!) — 2a alte lie! 
k>0 
Next we transform these densities to coordinate space 
and find for the local part 


(C1) 


p(rp) =2 DD 5 5 E E (minty Pain’ (PSP) 


m>0 nk s 


=> prlr.p) Pro 
Xr 


Using Eq. (A7) we finally obtain the various local densi- 
ties 


(C2) 


A= Eo nn (Gnu (Gua?) — Fan (t) F'n (r)] (em]Pa]k'm), (C3) 
mol nk 

E Dea 7 

RO) =D 5 A iat) [Gnn (r)Gnrw (r) + Fan(r)En w (r)] (sml Palk'm), (C4) 
m>0 nk 


and similar equations for the isovector density p3(r) and 
for the charge density p$ (r). 


Appendix D: Solution of the Klein-Gordon equation 


The various densities are the sources of the meson fields 
in the Klein Gordon equations (8). These equations are 
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solved by integrating the densities over the static Green 
functions in spherical coordinates. For simplicity we give 
here the details only for the o meson 


D(r,0,r',0, Mmo) = -Mo Y ja(imar<)hW (imers) 
X 
x (2A + 1) P; (cos 0) Py (cos 6’), (D1) 


and the photon 


A 
D(r,0,r',0') = Y -$ P, (cos 0) Py (cos 6’). 
Ae 


(D2) 


Here ry = max(r,r’) and re = min(r,r’). The solution 
for the ø field is 


= 5 oy(r)P (cos 8), 


A 


(D3) 


with 
oalr) = 4T JoMo (ratimer) f dr'jx limor’) p$ (r") 
0 


+ jalimor) i d'ha limor 05") , (D4) 


where j, and h are the spherical Bessel and Hankel func- 
tions. Similarly we find for the Coulomb field 


1 g 1 r'à pel AÀ Sa 1 1 Efa 
pt Jo E r T FINI py(7"). 
T: 


(D5) 
From the \ components of the meson fields ox (r), w9 (r) 
p9 (r) and AÌ (r) we find immediately the corresponding 
components of the scalar and the vector potential given 
in Eqs. (6) and (7). 


AX (r) = 


Appendix E: Pairing fields and tensors 


As in the case of the normal density we first calculate 
the pairing tensor « in the Woods Saxon basis 


(m)* m) 
= UE (mK) Uy, s(n’ Kn)" 


Next we transform it to coordinate space and obtain 

K(rsp,r’s’p'). This is a 2x2 matrix in spin space and 
therefore it can be expressed in terms of the unity and 
the Pauli matrices 


(E1) 


oE KJ 


k(rsp,r's'p') = K(rp,r'p') + r(rp,r'p')-o, (E2) 
where K(rp,r'p') is the S = 0 part and K(rp,r’p’) is a 
vector, the S = 1 part of the pairing tensor. We re- 
alize that the special form of the pairing interaction in 
Eq. (11) guarantees that we do not need the full ma- 
trix K(rsp,r's'p'). As a consequence of the zero range 
we need only the local part of this matrix, and since the 
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force acts in the S = 0 channel, only the spin scalar part 
of k contributes. It is obtained by coupling to S = 0: 


= J (-)+tk(rsp, r-sp'). 


S 


k(r, p,p") (E3) 


As mentioned above, in this work we take into account 
only pairing fields which are diagonal in the quantum 
number p. 

Because of the symplectic structure of the RHB equa- 
tions the pairing tensor «k connects basis states |nkm) 
with the time reversal states |n’«’m). Using the same 
arguments as in Eq. (B6) we obtain for the local and 
scalar part of the pairing density 


K(rp) = 2 5 > 5 op (rsp)k 


m>0 nk s 


=X (r, p)PA (O) 
A 


Finally with the help of Eq. (A7) we obtain the pairing 
densities in various -channels 


92 +1 
++ >D a 
0) En Gnn® (nn), (n ra) Gw K' 


m>0 nk 


N (rsp) 


(E4) 


x («Km|Py|K’m), (E5) 


ae, 52 +1 
Mya ee a > Fane), (n's) Enix! 


m>0 nk 


x (km|Pi|k'm). (E6) 


As a consequence of these simplifications the gap equa- 
tion (10) has the local form 


Vof (r)e(r, p) 


with f(r) = (1 — p(r)/psat). The decomposition of this 
equation into spherical harmonics yields 


A(r,p) = (E7) 


Py 2 
Alr) = (2\+1)Vo X` far(r yet) (99) 
a aA 
(E8) 


and f(r) = (dx0 — p(T) a/ Psat): 


Appendix F: Relativistic structure of the pairing 
field 


So far we have neglected parts of the pairing field which 
connect large and small components, i.e., we have as- 
sumed that 

Atv (r) =K (r)=0. (F1) 
Since the density p(r) and the density-dependent func- 
tion f(p(r)) of the pairing interaction does not mix these 
components, the structure of Eq. (E8) shows also very 
clearly that the pairing tensor «(r) and the pairing field 


201708.00252v1 


chinaXiv 


A(r) are in this respect completely connected. If «(r) 
mixes these components, so does A. 

Considering the structure of Eqs. (E5) and (E6) we 
find that a non-vanishing term x+- (r) would have the 
form 


92A +1 
+- 
PA (r) 2- ae) » 5 eTA A (n’k’) Fri! 


m>0 nk 


x (Km|Py|K'm).(F2) 


Since large and small components have different parity, 
non-vanishing values of this function would require odd 
values of | +l’ and because of the parity selection rule 
in (Km|P,|K'm) also odd values of A. This means the 
parts of kt~(r) and At~(r) can only be expanded in 
components with odd A values. Of course, this fact is 
rather trivial. It does not violate parity, because even 
the simple Dirac equation with parity conserving fields 
have large and small components with different parity: 
lin. =o-p has also a negative parity. 

We can conclude that in the spherical case, where 
à = 0 and even, the field At~(r) has to vanish. In 
the deformed case this is not necessarily true. On the 
other side, these considerations depend on the interac- 
tion, as for instance on the fact that the pairing force we 
have used here excludes S = 1. In particular we did not 
take into account odd A-values in the pairing field and 
therefore At~ fields are excluded from the beginning. If 
we would allow for S = 1 pairs, spin-vector components 
of the form A-a@ mixing large and small components 
are not excluded, even in the spherical case, because in 
this case L = 1 and S = 1 can couple to J = 0. Of 
course this depends on the interaction. In Ref. [82] the 
S = 1 part of the zero range pairing force was not ex- 
cluded and non-vanishing pairing fields A*~ were taken 
into account. However, they turned out to be an order 
of magnitude smaller than the diagonal matrix elements 
Att. In particular they are very small as compared to 
the term ø -p which mixes large and small components in 
the Dirac Hamiltonian. Therefore they can be neglected 
as a very good approximation. 


Appendix G: Microscopic center of mass correction 


The center of mass correction in Eq. (27) which is 
widely used in the literature [75] can be derived as a first 
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order correction to a projection after variation [108] onto 
good linear momentum. In Ref. [109] this term has been 
derived in the framework of the BCS approximation as 


h2 
Dm a O a o 


i>0 


+ 5 [vivi (Vivi + uiu) (Vil? + IVl] 
i,i'>0 

In the following we show, how this formula can be ap- 
plied in the framework of relativistic Hartree-Bogoliubov 
theory. In a first step we use the fact that any Hartree- 
Bogoliubov wavefunction can be expressed in the form 
of a BCS-state in the canonical basis [75]. This basis is 
obtained by the diagonalization of the density matrix in 
p = V*VT in the Woods-Saxon basis (18) 


m 2% 
` Pnr minib k! T Ui Cnr: 


n' K! 


(G2) 


The eigenvalues v? are the BCS-occupation probabilities 
and the eigenvectors are the expansion coefficients of the 
canonical wave functions in the spherical spinors of the 
Woods-Saxon basis 


(G3) 


= 5 CrunPnnm (18). 
NK 


Here, i = (nmr) where m is the third component of the 
angular momentum j and m = + is the parity. 


Of course the eigenvalues of Eq. (G2) provide us only 
the absolute values of the occupation amplitudes v; and 
ui = \/1— v?. In Eq. (G2) we also need the sign of uivi. 
It is determined by the diagonal elements of the pairing 
tensor in the canonical basis, 


Pu 


=) j bo at 
UjUi = CnkCn'K' *(nk)(n’k’) + 


NK 


(G4) 


For the direct term we need the diagonal matrix ele- 
ments of the Laplacian 


T> i Gnr (r) OGrin(r) , 
Au = 2 on ch f ar | EED PEt y qy p Eet 


(G5) 


For the exchange term, we have |V; |? = 5, (—)" (Vp )w (V—,)iv and according to the Wigner-Eckart theorem [107] 
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we find for the spherical coordinate u of the gradient operator 


Jw = ye Cutie ae 1 i) (dna l|V||Onin’)s (G6) 


mum 
with the reduced matrix element 


(nal rw) = (YF [ITT i Y fran MAT Gwt) 
HITTITE S AEO], (c7) 


where j = v2j +1 and the expression (I||V]|{l’) is the reduced matrix element of V with respect to the integration 
over the angles. Of course, it still contains derivatives with respect to the radial coordinate. Following Sect. 5.7 of 
Ref. [107], we obtain 


d l d V 
UTI = Sv vI Z = L — ôv- VU E + S : (G8) 
rT r T 


r 


= 


